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The CP-PACS and JLQCD Collaborations have been made systematic studies of lattice QCD carrying out both 
chiral and continuum extrapolations. Importance of dynamical quark effects has been clarified by a comparison of 
quenched QCD and two flavor full QCD simulations. In two flavor simulations, the dynamical effects of u and d 
quarks are taken into account, but the third quark s is still treated in a quenched approximation. As the final step 
towards a fully realistic lattice simulation of QCD, two collaborations started a joint project of three flavor QCD, 
concentrating all big computers available. Based on a series of preparative studies of exact simulation algorithm 
and non-perturbative improvement coefficient, a large scale simulation of three flavor QCD has been started. I 
present the results of light hadron spectrum and light quark masses from the first production runs. 



1. INTRODUCTION 

A direct calculation of hadronic properties 
from QCD is a fundamental objective in parti- 
cle physics. Strong correlations among quarks 
make a reliable analytic calculation difficult even 
for basic properties of hadrons such as the mass 
spectrum and decay constants. So far numeri- 
cal simulation based on the lattice formulation of 
QCD is the only reliable way toward this goal. 

A systematic study of lattice QCD carrying 
out both the chiral and continuum extrapola- 
tions requires a huge amount of computations, 
however. The CP-PACS and JLQCD Collabora- 
tions have been performing a series of large scale 
simulations of QCD adopting the fastest comput- 
ers available. The CP-PACS Collaboration have 
been mainly using the CP-PACS computer, a 
dedicated parallel computer with the peak perfor- 
mance of 614 GFLOPS. The CP-PACS was de- 
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signed and developed in 1996 at the Center for 
Computational Physics, University of Tsukuba 
PP. The main engine for the JLQCD Collabo- 
ration is the supercomputers at the High Energy 
Accelerator Research Organization (KEK). The 
central computer since 2000 is the Hitachi SR- 
8000/F1 with 100 nodes achieving the peak speed 
of 1.2 TFLOPS 0. 

As the first step, the CP-PACS Collaboration 
carried out an extensive study of QCD in the 
quenched approximation, in which the effects of 
dynamical quark pair creations and annihilations 
are suppressed Performing the first well- 

controlled chiral and continuum extrapolations 
on lattices with the spatial extent of about 3 fm, 
the limitation of the quenched approximation was 
clearly proven: the light hadron mass spectrum 
deviates from the experiment by O(10%). 

The next natural step is to incorporate the 
effects of dynamical u, d quarks (two flavor 
full QCD). The results of systematic studies by 
the CP-PACS and JLQCD Collaborations @E] 
show that the discrepancies in the hadron mass 
spectrum observed in the quenched study are 
largely reduced by the dynamical u, d quarks. 
This demonstrates the importance of dynamical 
quarks. 
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In two flavor simulations, however, the third s 
quark is treated in a quenched approximation yet. 
On the other hand, we do expect that contribu- 
tion of dynamical s quark is not small because its 
mass is of the same order of magnitude as a typi- 
cal energy scale of gluon dynamics: m s < Aqcd- 
Therefore, introduction of the dynamical s quark 
in the simulation is the last major step left to- 
wards a fully realistic simulation of QCD. 

The CP-PACS and JLQCD Collaborations 
have started a systematic study of three flavor 
full QCD. The first target is to perform a precise 
measurement of the light hadron spectrum and 
light quark masses. 

Because the required computer power is enor- 
mous, we have decided to concentrate all big com- 
puters available to us for this joint project. In Ta- 
ble ^ I list the major computers we are devoting 
to this project. Fractions of the peak performance 
which we may fully use for lattice QCD simula- 
tions are given for each computer. Summing up 
these numbers, we may use about 2.5 TFLOPS 
to simulate three flavor QCD. 

In this paper, I present the outline and sta- 
tus of the project. The first step was to de- 
velop and test an exact algorithm for three-flavor 
QCD simulations. From a comparison of vari- 
ous variants of exact algorithms as described in 
Sec. El we have chosen a polynomial hybrid Monte 
Carlo (PHMC) algorithm. In Sec. 01 I discuss our 
test study of three flavor QCD using the PHMC 
algorithm, which indicates that improvement of 
lattice action is essential for meaningful simula- 
tions. We thus adopt the RG-improved gauge ac- 
tion by Iwasaki and the clover-improved Wilson 
quark action. To achieve full 0(a)-improvement, 
we then determined the non-perturbative value of 
the clover coefficient csw- Based on these prepar- 
ative studies, we have now started a large scale 
simulation of three flavor QCD. Results of the 
first production runs are presented in Sec.0] Ten- 
tative conclusions are given in Sec. 

2. EXACT ALGORITHM 

With the Wilson-type (staggered-type) lattice 
quarks, the exact HMC algorithm exists only for 
the cases of even (four times integer) number of 
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Figure 1. Thermal cycles for plaquette obtained 
on 8 3 x 16 lattices at a" 1 ~ 1.5-2 GeV \TQ. a) 
Plaquette gauge action at (3 = 4.6, 4.8, ■ • ■, 6.0 
from bottom to top. b) RG gauge action at {3 — 
1.5-2.25 in steps of 0.05 from bottom to top. 



flavors. Because an exact algorithm for odd num- 
ber of flavors was not known before, inexact R- 
algorithm has been adopted in previous studies of 
three flavor QCD. To avoid possible unexpected 
systematic errors, however, it is highly desirable 
to adopt an exact algorithm for large scale simu- 
lations. 

Recently, several exact algorithms for odd fla- 
vors have been proposed: the multi-boson algo- 
rithm jS| and the polynomial hybrid Monte Carlo 
(PHMC) algorithm PTU| . both based on polyno- 
mial approximation for the quark matrix and for 
its determinant. We adopt the PHMC algorithm 
extending it to clover-improved Wilson quarks. 
We made a systematic test of several variants of 
the PHMC algorithm [TT]. We found that, with 
appropriate improvements, the PHMC algorithm 
is equally efficient as the HMC algorithm for the 
case of two-flavor QCD. 
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Table 1 



Computers for the three flavor QCD project by CP-PACS/JLQCD Collaborations. Performance of our 
PHMC code was measured in an actual production runs of three-flavor QCD on 20 3 x 40 lattices. 



machine 


location 


# nodes 


peak speed 
[GFLOPS] 


fraction for 
LQCD [GFLOPS] 


performance 
of PHMC code 


CP-PACS [J 


CCP, U.Tsukuba 


2048 


614 


-614 


20% 


SR-8000/G1 


CCP, U.Tsukuba 


12 


173 


-173 


44% 


SR-8000/F1 
VPP-5000 


KEK 


100 


1200 


-768 


35% 


SICP, U.Tsukuba 


80 


768 


-230 


44% 


Earth Simulator jj] 


ES Center 


640 


40960 


-640 


31% 



From this study, we decide to adopt the PHMC 
algorithm for the s quark and the HMC algo- 
rithm for degenerate u, d quarks. Note that, be- 
cause the s quark is heavy, the problem of neg- 
ative determinant at very small quark masses is 
automatically avoided. Together with additional 
improvements, we find that the overall CPU time 
to simulate one trajectory of three flavor QCD 
is only about 1.5 times more than that for two 
flavor QCD using our previous HMC code at the 
same u, d quark masses. We conclude that the 
PHMC algorithm is sufficiently efficient to carry 
out a systematic simulation of three flavor QCD 
with the present power of computers. 

We have implemented the code to the comput- 
ers listed in Tabled Optimizing the vectorization 
and parallelization algorithms depending on the 
characteristics of each machine, we achieved the 
performance of 20-44% in actual production runs 
[L2"] , as compiled in the last column of Table H 

3. CHOICE OF THE LATTICE ACTION 

Improvement of the lattice action is effective 
to suppress lattice artifacts on coarse lattices [T3"| 
and has played an essential role in our study of 
two flavor QCD to reduce necessary computer re- 
sources. Because the requirement of computer 
power is even more oppressive in the three fla- 
vor project due to larger number of parameters, 
improvement will be important for a systematic 
simulation of three flavor QCD. 

3.1. A test study 

Using the exact PHMC algorithm, we made 
a series of test studies at a -1 ~ 1.5-2 GeV on 



8 3 x 16 and 12 3 x 32 lattices [H]. For ghi- 
ons, we test the standard one-plaquette action, 
the RG- improved action by Iwasaki JHjj and 
the meanfield-improved Symanzik gauge action. 
For quarks, we adopt the clover-improved Wil- 
son quark action with meanfield-improved clover- 
coefficient csw- 

With the plaquette gauge action, we have en- 
countered a severe lattice artifact at a -1 < 2 
GeV: the plaquette expectation values shown in 
Fig. Ola) indicate unexpected first-order transi- 
tions at /3 = 4.95 and 5.0. Our study of the lat- 
tice size dependence suggest that this is a bulk 
transition. On the other hand, results with im- 
proved gauge actions show no signs of hysteresis, 
as shown in Fig. IHb) for the RG-improved ac- 
tion. We suspect that the lattice artifact in the 
case of plaquette action is due to an effective ad- 
joint coupling from the clover term |14| . 

In any case, because it is difficult to simulate 
several points beyond a -1 —2 GeV, our findings 
imply that improvement of the gauge action is 
indispensable to perform a continuum extrapola- 
tion. We adopt the RG-improved gauge action 
for gluons. 

3.2. Non-perturbative Csw 

For quarks, we adopt clover-improved Wilson 
quark action. To completely remove 0{a) errors, 
we need non-perturbative values of the clover co- 
efficient csw , which have not been estimated for 
three flavor QCD. In |16I17| . we have determined 
the non-perturbative csw both for the plaque- 
tte and RG-improved gauge actions. We use the 
Schrodinger functional method |18|. For the RG 
gauge action we adopt the boundary condition 
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Figure 2. Non-perturbative value of csw in three- 
flavor QCD with RG-improved gauge action ^7] . 
Open and filled triangles show csw at a fixed 
physical lattice size L* , obtained by correcting 
the raw results shown by open circles. Dotted 
line shows the one- loop result for L/a = oo. 



described in [T5]. 

We found that the finite volume effect in csw 
is not negligible for the case of the RG-improved 
gauge action. When we estimate csw at a fixed 
dimension- less lattice size L/a as usually done in 
previous studies, csw suffers a constant devia- 
tion from the true non-perturbative value even in 
the continuum limit. This yields 0(a) errors in 
physical observables. 

In order to remove the finite volume effect, 
we estimate the non-perturbative csw at a fixed 
physical lattice size L*. Then the finite volume 
effect vanishes as a/L* in the continuum limit. 
Correction of the data to L* from the results ob- 
tained at the simulation point L is done using 
one-loop formulae calculated with the SF setup 
in finite volume [HJ\. See J7| for details. 

Our final results for csw at L* = 6a^ =6 / s 2 =1 9 
are summarized in Fig. [2] by triangles. We note 
that the results smoothly converge to the one- 
loop result for L/a = oo (dashed line) in the con- 
tinuum limit go = 0. 

4. THREE FLAVOR QCD 

Based on the preparative studies discussed in 
previous sections, we have started a systematic 
simulation of three flavor QCD. As the first point 
towards the continuum limit, we are carrying out 
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Figure 3. Simulation points on 16 3 x 32 lattice at 
(3 = 1.9. Our estimation for the physical point is 
given by the star. The line shows the chiral limit 
defined by mpg,LL = 0. 

jobs at a -1 ~ 2 GeV (a ~ 0.1 fm). We have fin- 
ished the production runs on a 16 3 x 32 lattice 
(La~ 1.6 fm), and are now carrying out simula- 
tions on a 20 3 x 40 lattice (La ^2.0 fm). In this 
section, I present the results of the first produc- 
tion runs on the 16 3 x 32 lattice [201 - 

4.1. Simulation parameters 

With the RG-improved gauge action and 
clover-improved Wilson quark action, we made 
simulations at (3 — 1.9 on a 16 3 x 32 lattice |2T]|. 
The non-perturbative value of csw at this (3 is 
1.715 (TJ|. We studied six values of the u, d quark 
mass in the range K u d — 0.1358-0.1370 corre- 
sponding to TOps,LL/f7iv,LL — 0.64-0.77, where 
L means the light u, d sea quark. For the s 
quark mass, we studied two points K s = 0.1364 
and 0.1358 corresponding to mps,ss/wv,ss — 0.72 
and 0.77, where S is for the sea s quark. These 
values are close to the physical s quark point 
to^ s /to0~O.68 from the chiral perturbation the- 
ory. The simulation points are summarized in 
Fig. El 

We simulated 3000 trajectories at each 
(Kud, K s ) and accumulated the configurations ev- 
ery 10 trajectories. The HMC step size and the 
order of polynomial in the PHMC algorithm were 
adjusted to achieve the acceptance rate of more 
than 85% in HMC steps and more than 90% in 
PHMC steps. So far, we have measured hadronic 
observables only when each valence quark is one 
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Figure 4. Vector meson mass tov.ll at (3 = 1.9 
on 16 3 x 32 lattice. 



of the sea quarks, L or S. Errors are estimated by 
a jack-knife method with bins of 50 trajectories. 

4.2. Meson spectrum 

Because the spatial lattice size of La ~ 1.6 fm 
is not quite large, we concentrate on mesons in 
this report. We tested point and exponentially 
smeared sources. Because we obtained clearest 
plateaus for effective masses when both quark 
sources are smeared, we show the results with 
doubly smeared sources in the following. 

Results for the mass of vector mesons consist- 
ing of two light sea quarks are shown in Fig. 0] 
as a function of 1/K — 1/K C , where K c is de- 
fined by mps,Lh(K u d — K s = K c ) = 0. Results for 
other mesons are similar. Because the quark mass 
dependence is quite smooth, we adopt following 
polynomial ansatze 
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Figure 5. 4> meson mass from the if-input as a 
function of the lattice spacing. The result from 
three flavor QCD is given by the filled circle at 
a~0.1 fm. Quenched (open symbols) |3I4| and 
two flavor QCD results (shaded symbols) |4I5| are 
also shown. 



both combinations of JU+GJ and 0+©- Both 
ansatze fit the data well. We quote the weighted 
average of two combinations as the central value 
of the masses, while the difference between two 
combinations is treated as a systematic error from 
the chiral fit. To identify the lattice spacing a 
and the physical point (K^ ys , Kf hvs ), we use ei- 
ther (M„,M p ,M K ) (X-input) or (M n ,M p ,M<f,) 
(0-input). 

Result for from the -ftf-input is shown in 
Fig-El The clear discrepancy between experiment 
(star) and quenched results (open symbols) is 
largely removed in two flavor QCD (shaded sym- 
bols) in the continuum limit. Our new three fla- 
vor value is shown by the filled symbol at a ~ 0.1 
fm. Wc find that the three flavor result is higher 
than the two flavor results and, already at a ~ 0.1 
fm, the mass is consistent with experiment. Re- 
sult for Mk* given in |2()| is similar. 

We summarize the relative discrepancy with 
experiment in Fig. [B| We see that experimental 
spectrum is well reproduced in three flavor QCD 
at a ~ 0.1 fm. Reproduction of experimental val- 
ues can be confirmed also by the J parameter 
20 for which a chiral fit is not required. Consis- 
tency with experiment at finite a imply that the 
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Figure 6. Relative discrepancy of meson masses 
in three flavor QCD at a" 1 ~2 GeV 20 . Results 
of quenched QCD in the continuum limit are also 
shown by open symbols. 



K- and 0-inputs lead to consistent values already 
at finite a. Accordingly, we obtain a -1 = 2.05(4) 
GeV from the if- input and 2.05(5) GeV from the 
0- input. 

For a precise prediction, we need to extrapo- 
late the results to the continuum limit. However, 
we may hope that the scaling violation is indeed 
small for our choice of non-perturbatively 0(a)- 
improved action, such that the consistency with 
experiment is kept in the continuum limit. This 
should be tested in future. 

4.3. Quark mass 

Because quarks are confined, their masses are 
not direct observables of the theory. Therefore, 
there exist several alternative definitions for the 
quark mass. Two popular definitions are the 
axial-vector Ward identity (AWI) quark mass 
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(4) 



and the vector Ward identity (VWI) quark mass 



m q = Z'(l/K-l/K c )/2, 



(5) 



where A 4 is the fourth component of the axial- 
vector current, P is the pseudo-scalar density, and 
Z's are renormalization factors. Different defini- 
tions lead to different values for m q on finite lat- 
tices. This was a big source of error in early calcu- 
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Figure 7. Quark masses (a) m u d and (b) m s in 
the MS scheme at n = 2 GeV |2D] . 



lations of m q . In our previous studies of quenched 
and two flavor full QCD, however, we have shown 
that they converge to universal values in the con- 
tinuum limit |3I4| . 

In three flavor QCD, we found that, although 
the differences between K- and </>-inputs are ab- 
sent, the AWI and VWI quark masses disagree at 
a -1 ~ 2 GeV. We also noted that the choice of 
"ifc" in © introduces a sizable ambiguity in the 
values of VWI quark mass: From Fig. we ob- 
tain even negative value of m u d when we adopt K c 
defined by mps t hh(K u d = K S = K C ) = 0. Such am- 
biguities will be removed in the continuum limit. 

At the present stage having data only at a value 
of a -1 , we would like to concentrate on the AWI 
quark mass which does not suffer from the ambi- 
guity of K c and, in the case of two flavor QCD, 
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shows a smaller scaling violation than VWI quark 
masses 0]. To convert to the quark mass in the 
MS scheme, we match the renormalized mass with 
the lattice data at [i = a^ 1 using a meanfield- 
improved one-loop Z factor |21j , and let it run to 
/i = 2 GeV using the 4-loop beta function. 

Our results are summarized in Fig.0 Previous 
results from quenched and two flavor full QCD, 
extrapolated to the continuum limit, are also 
shown for comparison |3I4| . We see that quark 
masses decrease when we increase the number of 
dynamical quark flavors. In three flavor QCD, we 
obtain m^f = 2.89(6) MeV and mf^ = 75.6(3.4) 
MeV at a -1 ~2 GeV, where the central values are 
from the if-input and the errors include system- 
atic errors estimated by the difference between 
the K- and 0-inputs. With the effects of the dy- 
namical s quark, both m ue i and m s are decreased 
by about 15% from the previous two flavor val- 
ues. The ratio m s /m uc i — 26.2(1.0) is consistent 
with the one-loop estimate of chiral perturbation 
theory 24.4(1.5) E2| 



5. TENTATIVE CONCLUSIONS AND 
OUTLOOK 

We have presented the status of the joint 
project of three flavor QCD simulation by the 
CP-PACS and JLQCD Collaborations. We em- 
ploy exact HMC algorithm for light u, d quarks 
and exact PHMC algorithm for s quark. To 
reduce the lattice artifacts, we adopt the RG- 
improved gauge action and non-perturbatively 
0(a)-improved clover quark action. As the first 
step toward a systematic study of the fully real- 
istic QCD on the lattice, we made simulations at 
a" 1 ~ 2 GeV on 16 3 x 32 lattices. These sim- 
ulations show that, with the effects of dynami- 
cal s quark, light meson mass spectrum agrees 
well with experiment already at a -1 ~ 2 GeV. 
We also found that the quark masses m u d and 
m s are lower than the two flavor values by about 
15%. We are currently performing simulations at 
the same simulation points on 20 3 x 40 lattices 
to study finite lattice volume effects. Although a 
continuum extrapolation is not made yet, these 
results are quite encouraging to further carry out 



large scale simulations of three flavor QCD. 
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